GetXYSizesFromFile Subroutine

private subroutine GetXYSizesFromFile(ncId, x, y, idx, idy, shape, latlon)

extraxt the number of columns (x) and rows (y) from netCDF file.

Arguments

Type IntentOptional Attributes Name
integer(kind=short), intent(in) :: ncId

NetCdf Id for the file
integer(kind=short), intent(out) :: x

number of columns and rows
integer(kind=short), intent(out) :: y

number of columns and rows
integer(kind=short), intent(out), optional :: idx

id of x variable
integer(kind=short), intent(out), optional :: idy

id of y variable
character(len=*), intent(out) :: shape

ccordinate in vector or matrix
integer(kind=short), intent(out), optional :: latlon

coordinate order is lat-lon (1) or lon-lat (2)

Variables

Type Visibility Attributes Name Initial
integer(kind=short), public :: IdXCoordinate

Id of the variable containing information on x ccordinate
integer(kind=short), public :: IdXDim

Id of the X dimension

integer(kind=short), public :: IdYCoordinate

Id of the variable containing information on y ccordinate
integer(kind=short), public :: IdYDim

Id of the Y dimension

character(len=80), public :: attribute
integer(kind=short), public, ALLOCATABLE :: dimIds(:)

dimension Ids of a variable
integer(kind=short), public :: i
logical, public :: longnamefound
integer(kind=short), public :: nDimsVar

number of dimensions of a variable
integer(kind=short), public :: nVars

number of variables
integer(kind=short), public :: ncStatus

error code return by NetCDF routines
logical, public :: stdnamefound
logical, public :: unitsfound

Source Code

SUBROUTINE GetXYSizesFromFile &
!
(ncId, x, y, idx, idy, shape, latlon)

IMPLICIT NONE

!Arguments with intent(in):
INTEGER (KIND = short), INTENT(IN) :: ncId  !!NetCdf Id for the file

!Arguments with intent(out):
INTEGER (KIND = short), INTENT(OUT) :: x, y !!number of columns and rows
INTEGER (KIND = short), OPTIONAL, INTENT(OUT) :: idx !! id of x variable
INTEGER (KIND = short), OPTIONAL, INTENT(OUT) :: idy !! id of y variable
CHARACTER (LEN = *), INTENT(OUT) :: shape !!ccordinate in vector or matrix
INTEGER (KIND = short), OPTIONAL, INTENT(OUT) :: latlon !! coordinate order is lat-lon (1) or lon-lat (2)


! Local variables:
INTEGER (KIND = short)          :: ncStatus !!error code return by NetCDF routines
INTEGER (KIND = short)          :: nVars !!number of variables
INTEGER (KIND = short)          :: IdXCoordinate !!Id of the variable containing 
                                                 !!information on x ccordinate
INTEGER (KIND = short)          :: IdYCoordinate !!Id of the variable containing 
                                                 !!information on y ccordinate
INTEGER (KIND = short)          :: IdXDim !!Id of the X dimension 
INTEGER (KIND = short)          :: IdYDim !!Id of the Y dimension 
INTEGER (KIND = short)          :: nDimsVar !!number of dimensions of a variable
INTEGER (KIND = short), ALLOCATABLE :: dimIds (:) !!dimension Ids of a variable
CHARACTER (LEN = 80)            :: attribute
INTEGER (KIND = short)          :: i
LOGICAL                         :: stdnamefound
LOGICAL                         :: longnamefound
LOGICAL                         :: unitsfound

!------------end of declaration------------------------------------------------


!inquire dataset to retrieve number of dimensions, variables 
!and global attributes
ncStatus = nf90_inquire(ncId, nVariables = nVars )
                  
CALL ncErrorHandler (ncStatus)

stdnamefound = .FALSE.
longnamefound = .FALSE.
unitsfound = .FALSE.

!scan dataset searching for variable referring to grid mapping
DO i = 1, nVars
  attribute = ''
  ncStatus = nf90_get_att (ncId, varid = i, name = 'standard_name', &
                           values = attribute)
  
  IF (ncStatus == nf90_noerr) THEN !'standard_name' is  found
    IF ( attribute(1:23) == 'projection_x_coordinate' .OR. & !projected reference system
         attribute(1:9) == 'longitude' .OR. & !geographic reference system
         attribute(1:7) == 'Easting' .OR. & 
         attribute(1:14) == 'grid_longitude')THEN !rotated pole system
       stdnamefound = .TRUE.
       IdXCoordinate = i
         ELSE IF ( attribute(1:23) == 'projection_y_coordinate' .OR. &
              attribute(1:8) == 'latitude' .OR. & !geographic reference system
              attribute(1:8) == 'Northing' .OR. & 
              attribute(1:13) == 'grid_latitude')THEN !rotated pole system
       stdnamefound = .TRUE.
       IdYCoordinate = i
    END IF
  END IF

END DO

IF ( .NOT. stdnamefound) THEN
    CALL Catch ('warning', 'GridLib',        &
    'standard name not found while searching for x and y axis')
END IF

DO i = 1, nVars
  attribute = ''
  
  ncStatus = nf90_get_att (ncId, varid = i, name = 'long_name', &
                           values = attribute)

  IF (ncStatus == nf90_noerr) THEN !'long_name' is  found
    IF ( attribute(1:23) == 'projection_x_coordinate' .OR. & !projected reference system
         attribute(1:9) == 'longitude' .OR. & !geographic reference system
         attribute(1:7) == 'Easting' .OR. & 
         attribute(1:14) == 'grid_longitude')THEN !rotated pole system
       longnamefound = .TRUE.
       IdXCoordinate = i
         ELSE IF ( attribute(1:23) == 'projection_y_coordinate' .OR. &
              attribute(1:8) == 'latitude' .OR. & !geographic reference system
              attribute(1:8) == 'Northing' .OR. &
              attribute(1:13) == 'grid_latitude')THEN !rotated pole system
       longnamefound = .TRUE.
       IdYCoordinate = i
    END IF
  END IF

END DO

IF ( .NOT. longnamefound) THEN
    CALL Catch ('warning', 'GridLib',        &
    'long_name not found while searching for x and y axis')
END IF

DO i = 1, nVars
  attribute = ''

  ncStatus = nf90_get_att (ncId, varid = i, name = 'units', &
                           values = attribute)

  IF (ncStatus == nf90_noerr) THEN !'units' is  found
    IF ( attribute(1:11) == 'degree_east' )THEN 
       unitsfound = .TRUE.
       IdXCoordinate = i
    ELSE IF ( attribute(1:12) == 'degree_north' )THEN 
       unitsfound = .TRUE.
       IdYCoordinate = i
    END IF
  END IF
 
END DO

IF ( .NOT. unitsfound) THEN
    CALL Catch ('warning', 'GridLib',        &
    'units not found while searching for x and y axis')
END IF

IF ( .NOT. stdnamefound .AND. &
     .NOT. longnamefound .AND. &
      .NOT. unitsfound ) THEN
     CALL Catch ('error', 'GridLib',        &
    'x and y axis can not be found')
END IF

IF (PRESENT (idx) ) THEN
  idx = IdXCoordinate
END IF

IF (PRESENT (idy) ) THEN
  idy = IdYCoordinate
END IF

!retrieve Id of X  dimension
ncStatus = nf90_inquire_variable (ncId, IdXCoordinate, ndims = nDimsVar)
CALL ncErrorHandler (ncStatus)
ALLOCATE (dimIds (nDimsVar))
ncStatus = nf90_inquire_variable (ncId, IdXCoordinate, dimids = dimIds)
CALL ncErrorHandler (ncStatus)
IF ( nDimsVar == 1) THEN !x is a vector
  IdXDim = dimIds (1)
  shape = 'v'
ELSE IF ( nDimsVar == 2) THEN 
  !x is a matrix. Used for not regular grids.
  !Suppose that x is the first dimension
  IdXDim = dimIds (1)
  shape = 'm'
END IF
DEALLOCATE (dimIds)

!retrieve Id of Y dimension
ncStatus = nf90_inquire_variable (ncId, IdYCoordinate, ndims = nDimsVar)
CALL ncErrorHandler (ncStatus)
ALLOCATE (dimIds (nDimsVar))
ncStatus = nf90_inquire_variable (ncId, IdYCoordinate, dimids = dimIds)
CALL ncErrorHandler (ncStatus)
IF ( nDimsVar == 1) THEN !y is a vector
  IdYDim = dimIds (1)
ELSE IF ( nDimsVar == 2) THEN 
  !y is a matrix. Used for not regular grids.
  !Suppose that y is the second dimension
  IdYDim = dimIds (2)
END IF
DEALLOCATE (dimIds)

!set coordinate order convention lon-lat or lat-lon
IF (PRESENT (latlon) ) THEN
    IF (IdXDim > IdYDim) THEN !coordinate convention order is lon-lat
        latlon = 1
        !latlon = 2  !DAO
    ELSE !coordinate convention order is lat-lon
        latlon = 2
        !latlon = 1  !DAO
    END IF
END IF

!retrieve length
ncStatus = nf90_inquire_dimension (ncId, IdXDim, len = x)
CALL ncErrorHandler (ncStatus)

ncStatus = nf90_inquire_dimension (ncId, IdYDim, len = y)
CALL ncErrorHandler (ncStatus)

RETURN
END SUBROUTINE GetXYSizesFromFile